1. （自分のPCまたは教室のPCに）ログイン

2. ウェッブ・ブラウザー（Google Chrome など）を起動

   • Moodle の GES001 経済と経済学のサイトから、このスライドのページを表示（リンク「Rでデータサイエンス」https://ds-sl.github.io/intro2r/ges001/）

3. （別のタブまたは ウィンドウで）PositCloud にログイン[Posit.cloud]

   • アカウントのない人はサイン・アップ [共有プロジェクト] から、Save a Permanent Copy）

   • RStudio を自分のコンピュータにインストールしている人は起動

4. リンクの右上の Raw ボタンの右の Copy a raw file からコピーして演習用 R Markdown ファイルを作成（あとで再度解説します）[Rmd]

01/18(TH)　南部アフリカ諸国の貧困と不平等に対する対策

　　　　　　COVID-19が貧困に与えた影響 　　　　

• 機会の不平等と成果の不平等が相互に関係し、長期にわたって持続している。

• 人生の初期に出現する不平等を最小化するための政策が重要である。

• 要因：教育、労働市場、失業率、男女間格差、土地所有権、などなど　　

01/23(TU)　Rでデータサイエンス5　 [Main]

• Poverty headcount ratio at $6.85 a day (2017 PPP) (% of population)：SI.POV.UMIC [Link]

• Government expenditure on education, total (% of GDP)：SE.XPD.TOTL.GD.ZS [Link]

• School enrollment, primary (% gross)：SE.PRM.ENRR [Link]

• School enrollment, secondary (% gross)：SE.SEC.ENRR [Link]

• School enrollment, tertiary (% gross)：SE.TER.ENRR [Link]

• Mortality rate, under-5 (per 1,000 live births)：SH.DYN.MORT [Link]

• Incidence of HIV (% of uninfected population ages 15-49)：SH.HIV.INCD.ZS [Link]

• School enrollment, primary and secondary (gross), gender parity index (GPI)：SE.ENR.PRSC.FM.ZS [Link]

• Ratio of female to male labor force participation rate (%) (modeled ILO estimate)：SL.TLF.CACT.FM.ZS [Link]

• Unemployment, female (% of female labor force) (modeled ILO estimate)：SL.UEM.TOTL.FE.ZS [Link]

• Unemployment, male (% of male labor force) (modeled ILO estimate)：SL.UEM.TOTL.MA.ZS [Link]

• Net official development assistance and official aid received (current US$) DT.ODA.ALLD.CD [Link]

準備

library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.4.4     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.0
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(WDI)

データの読み込み（importing）

df_poverty_inequality <- WDI(
  indicator = c(gini = "SI.POV.GINI",
                under_6.85 = "SI.POV.UMIC",
                ed_exp = "SE.XPD.TOTL.GD.ZS",
                primary = "SE.PRM.ENRR",
                secondary = "SE.SEC.ENRR",
                tertiary = "SE.TER.ENRR",
                under5 = "SH.DYN.MORT",
                new_hiv = "SH.HIV.INCD.ZS",
                school_gpi = "SE.ENR.PRSC.FM.ZS",
                job_gpi = "SL.TLF.CACT.FM.ZS",
                female_unemploy = "SL.UEM.TOTL.FE.ZS",
                male_unemploy = "SL.UEM.TOTL.FE.ZS",
                oda = "DT.ODA.ALLD.CD"), extra = TRUE)

保存と読み込み

２回目からは、data から読み込めるようにしておきます。

最初の１回目は、かならず実行してください。

write_csv(df_poverty_inequality, "data/poverty_inequality.csv")

df_poverty_inequality <- read_csv("data/poverty_inequality.csv")

## Rows: 16758 Columns: 24
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr   (7): country, iso2c, iso3c, region, capital, income, lending
## dbl  (15): year, gini, under_6.85, ed_exp, primary, secondary, tertiary, und...
## lgl   (1): status
## date  (1): lastupdated
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

データを見てみよう (viewing)

df_poverty_inequality 

変数の選択（selecting）

df_pov_ineq <- df_poverty_inequality |> 
  select(country, iso2c, year, gini:oda, region, income, lending)
df_pov_ineq

変数同士の相関関係: NA ではないところのみ選択

df_pov_ineq |> drop_na(gini:oda) |> select(gini:oda)

cor(cars$speed,cars$dist)

## [1] 0.8068949

cars |> ggplot(aes(speed, dist)) + geom_point() + 
  geom_smooth(formula = 'y~x',method = "lm", se=FALSE)

相関係数：直線の傾きが正なら正、負なら負、直線に近い程、1 または-1 に近い

変数相互の相関係数を一度に求める

df_pov_ineq |> drop_na(gini:oda) |> select(gini:oda) |> cor() |> 
  round(digits = 2) |> as.data.frame()

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(gini, primary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(gini,primary) = 0.34")

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(gini, tertiary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(gini,tertiary) = -0.29")

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(under_6.85, tertiary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(under_6.85,tertiary) = -0.77")

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(under_6.85, under5)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(under_6.85,under5) = 0.71")

以下の指標の中から、一つを選択して、データの概要（description）を記録し、データを WDI で取得し、以下の分析をする。

1. 各年毎のデータの数の棒グラフ
2. 日本のデータの年の降順での表示
3. 経年変化を表す折れ線グラフ
   1. 日本
   2. 南部アフリカ関税同盟の５カ国
   3. 選択したいくつかの国
4. データが十分ある最近の年の値のヒストグラム
5. データが十分ある最近の年の値の10カ国の値の棒グラフ
   1. 値が大きい方から
   2. 値が小さい方から

それぞれについて考察（気づいたこと、疑問など）を記す

2023.1.25. 23:59 までに Moodle の演習の課題ボックスに提出したものについては、なるべく、早く見て、フィードバックを書きます。それ以降に提出されたものも見ますが、フィードバックは遅くなると思ってください。

データ

1. Government expenditure on education, total (% of GDP)：SE.XPD.TOTL.GD.ZS [Link] 変数名：ed_exp

2. School enrollment, primary (% gross)：SE.PRM.ENRR [Link] 変数名：primary

3. School enrollment, secondary (% gross)：SE.SEC.ENRR [Link] 変数名：secondary

4. School enrollment, tertiary (% gross)：SE.TER.ENRR [Link] 変数名：tertiary

5. Mortality rate, under-5 (per 1,000 live births)：SH.DYN.MORT [Link] 変数名：under5

6. Incidence of HIV (% of uninfected population ages 15-49)：SH.HIV.INCD.ZS [Link] 変数名：new_hiv

7. School enrollment, primary and secondary (gross), gender parity index (GPI)：SE.ENR.PRSC.FM.ZS [Link] 変数名：school_gpi

8. Ratio of female to male labor force participation rate (%) (modeled ILO estimate)：SL.TLF.CACT.FM.ZS [Link] 変数名：job_gpi

9. Unemployment, female (% of female labor force) (modeled ILO estimate)：SL.UEM.TOTL.FE.ZS [Link] 変数名：female_unemploy

10. Unemployment, male (% of male labor force) (modeled ILO estimate)：SL.UEM.TOTL.MA.ZS [Link] 変数名：male_unemploy

11. Net official development assistance and official aid received (current US$) DT.ODA.ALLD.CD [Link] 変数名：oda

スライドで見ている方は、RNotebook ファイルで見てください [リンク]

概要：国内総生産（GDP）に対する、国の教育に関する支出（Government expenditure on education, total (% of GDP)）のデータの分析を行う

Government expenditure on education, total (% of GDP)：SE.XPD.TOTL.GD.ZS [Link]

データ情報

• データ名：

• データコード：

• 変数名：

• 概要：

データの取得

準備

library(tidyverse)
library(WDI)

WDI パッケージを使って、直接データをダウンロードし、変数名を、ed_exp に指定。

df_ed_exp <- WDI(indicator = c(ed_exp = "SE.XPD.TOTL.GD.ZS"))

write_csv(df_ed_exp, "data/ed_exp.csv")

df_ed_exp <- read_csv("data/ed_exp.csv")

## Rows: 16758 Columns: 5
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (3): country, iso2c, iso3c
## dbl (2): year, ed_exp
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

データの確認

df_ed_exp

str(df_ed_exp)

## spc_tbl_ [16,758 × 5] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
##  $ country: chr [1:16758] "Africa Eastern and Southern" "Africa Eastern and Southern" "Africa Eastern and Southern" "Africa Eastern and Southern" ...
##  $ iso2c  : chr [1:16758] "ZH" "ZH" "ZH" "ZH" ...
##  $ iso3c  : chr [1:16758] "AFE" "AFE" "AFE" "AFE" ...
##  $ year   : num [1:16758] 2022 2021 2020 2019 2018 ...
##  $ ed_exp : num [1:16758] 3.91 4.63 4.35 4.54 4.74 ...
##  - attr(*, "spec")=
##   .. cols(
##   ..   country = col_character(),
##   ..   iso2c = col_character(),
##   ..   iso3c = col_character(),
##   ..   year = col_double(),
##   ..   ed_exp = col_double()
##   .. )
##  - attr(*, "problems")=<externalptr>

国と地域

country には、国と地域両方が入っています。地域の iso2c は以下のものです。

REGION <- c("1A", "1W", "4E", "7E", "8S", "B8", "EU", "F1", "OE", "S1", 
"S2", "S3", "S4", "T2", "T3", "T4", "T5", "T6", "T7", "V1", "V2", 
"V3", "V4", "XC", "XD", "XE", "XF", "XG", "XH", "XI", "XJ", "XL", 
"XM", "XN", "XO", "XP", "XQ", "XT", "XU", "XY", "Z4", "Z7", "ZF", 
"ZG", "ZH", "ZI", "ZJ", "ZQ", "ZT")

df_ed_exp |> filter(iso2c %in% REGION) |> distinct(country, iso2c)

df_ed_exp |> filter(!(iso2c %in% REGION)) |> distinct(country, iso2c)

分析する国のリスト

南部アフリカ関税同盟 The Southern African Customs Union (SACU)

SOUTH_AFRICA_FIVE <- c("South Africa", "Namibia", "Eswatini", "Botswana", "Lesotho")

ラテンアメリカでジニ指数が大きい４カ国

CHOSEN_GINI_COUNTRIES <- c("Suriname", "Belize", "Brazil", "Colombia")

1. 各年毎のデータの数の棒グラフ

df_ed_exp |> drop_na(ed_exp) |> filter(!(iso2c %in% REGION)) |>
  ggplot(aes(year)) + geom_bar()

2. 日本の教育費（% of GDP）

df_ed_exp |> filter(country == "Japan") |> 
  drop_na(ed_exp) |> arrange(desc(year))

3. 経年変化

a. 日本

df_ed_exp |> filter(country == "Japan") |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line()

気づいたこと・疑問

• 1970年代の急激な上昇、1990年ごろの急激な現象は、何が原因なのだろう。

• 2014年ごろから減少、2018年ごろから増加、2020年から2021年は減少。

b. 南部アフリカ関税同盟

df_ed_exp |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country))

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

df_ed_exp |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

気づいたこと・疑問

• 平均で見ると、上昇してきており、7% 程度という大きな割合になっている。

c. ラテンアメリカ４カ国

df_ed_exp |> filter(country %in% CHOSEN_GINI_COUNTRIES) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country))

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

df_ed_exp |> filter(country %in% CHOSEN_GINI_COUNTRIES) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

分布

データの数から、まずは、2020年について見てみる。

df_ed_exp |> filter(year == 2020) |> filter(!(country %in% REGION))|>
  drop_na(ed_exp) |>
  ggplot(aes(ed_exp)) + geom_histogram(binwidth = 1)

参考：SACU の５カ国の値を縦線で書き込むには下のようにします。

df_ed_exp |> filter(year == 2020) |> filter(country %in% SOUTH_AFRICA_FIVE) 

参考：日本とSACU の５カ国の値を縦線で書き込むには下のようにします。

JP <- 3.416981
SAF <- df_ed_exp |> filter(year == 2020) |> filter(country %in% SOUTH_AFRICA_FIVE) |> pull(ed_exp)
df_ed_exp |> filter(year == 2020) |> filter(!(country %in% REGION))|>
  drop_na(ed_exp) |>
  ggplot() + geom_histogram(aes(ed_exp), binwidth = 1) +
  geom_vline(xintercept = SAF, col = "red") + geom_vline(xintercept = JP, col = "blue") +labs(title = "2020年の教育費の対GDP百分率", subtitle = "日本：青、SACU：赤")

データが十分ある最近の年の値の10カ国の値の棒グラフ

a. 値が大きい方から

df_ed_exp |> filter(year == 2020) |> drop_na(ed_exp) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(desc(ed_exp)) |> head(10) |> 
  ggplot(aes(fct_reorder(country, ed_exp), ed_exp)) + geom_col() + 
  coord_flip() + labs(title = "Top 10 Countries", x = "country", y = "Government expenditure on education, total (% of GDP)")

b. 値が小さい方から

df_ed_exp |> filter(year == 2020) |> drop_na(ed_exp) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(ed_exp) |> head(10) |> 
  ggplot(aes(fct_rev(fct_reorder(country, ed_exp)), ed_exp)) + geom_col() + 
  coord_flip() + labs(title = "Lowest 10 Countries", x = "country", y = "Government expenditure on education, total (% of GDP)")

演習の内容と課題

基本的には、PositCloud（https://posit.cloud/）を使って実習

• 探索的データ分析（EDA） -
  • 練習と一つ目の課題（w5eda.Rmd） [リンク], [Rmd]
  • 二つ目以降の課題：w5eda1.Rmd [リンク], [Rmd]
• 課題：2023.1.25. 23:59 までに Moodle の演習の課題ボックスに提出したものについては、なるべく、早く見て、フィードバックを書きます。それ以降に提出されたものも見ますが、フィードバックは遅くなると思ってください。

1. 「みんなのデータサイエンス - Data Science for All」[はじめてのデータサイエンス]

   • 導入として、GDP（国内総生産）のデータを使って説明しています。

2. Posit Recipes（旧 Posit Primers）: The Basics 対話型の演習サイトの最初 [Link]

3. Posit Cheat Sheet. 早見表です。印刷して使うために、PDF も提供しています。[Site Link]

4. DataCamp Cheat Sheet: Tidyverse for Biginners. データサイエンスの教育をしている会社の早見表の一つです。基本が簡単にまとまっています。[Link]